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Abstract 

Three regions of excessive flux of cosmic rays with energies of the order of PeV are found in 
the experimental data of the EAS MSU array at a confidence level greater than 4(T. For two of 
them, there are similar regions in the experimental data of the EAS-1000 Prototype array. One 
yS^ . of the interesting features of the regions is the absence of supernova remnants in their vicinities, 

Qh' traditionally considered as the main sources of Galactic cosmic rays, but the presence of isolated 

Q ■ pulsars, some of which are able to accelerate heavy nuclei up to energies close to PeV. In our opinion, 

H ' this favors the assumption that isolated pulsars are able to contribute to the flux of Galactic cosmic 

C/3 , rays more than is usually assumed. 

CN ■ 1 Introduction 

>, 

The problem of the origin of cosmic rays (CRs) in the PeV energy range remains open for several 
f-~. . decades. The model of diffusive shock acceleration at the outer front of expanding supernova remnants 

CN \ (SNRs) has become the most widely spread one, see, e.g., [1], though it has not been experimentally 

C^ ■ confirmed yet. Along with supernova remnants, other possible sources of PeV CRs have been con- 

sidered. In particular, shortly after pulsars were discovered and identified with neutron stars, it was 
demonstrated that they can be effective accelerators of CRs [2]. Interest to pulsars as possible sources 
of CRs with energies > 10^^ eV and up to ultrahigh energies did not vanish [3-7]. Models considered 
included acceleration both in pulsar wind nebulae (e.g., the Crab nebula) and by isolated pulsars. It 
was demonstrated in particular that the Geminga pulsar is a possible candidate for being the single 
C^ ' source of the knee at 3 PeV [5]. There is no surprise then that one of the areas of investigations in 

the field are attempts to find correlation between arrival directions of CRs and coordinates of their 
possible astrophysical sources, including SNRs, pulsars, and some other. 

In our previous works dedicated to the analysis of arrival directions of extensive air showers (EAS) 
registered with the EAS MSU and EAS-1000 Prototype (below, "PRO-1000") arrays, we have already 
presented a number of regions of excessive fiux (REFs) of PeV cosmic rays [8-11]. It was shown that 
not all of them could be matched with Galactic supernova remnants. However, a substantial part of 
REFs could be matched with isolated pulsars in whose vicinities nebulae are not observed. In this 
work, we consider three regions selected in the EAS MSU array data at a level of confidence greater 
than 4(7. None of these could be matched with Galactic supernova remnants, but there are isolated 
pulsars inside or in close vicinities of these regions, some of which are able to accelerate heavy nuclei 
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up to energies close to PeV. For two of these regions, there are similar regions in the PRO-1000 
experimental data. 

2 Experimental Data and Method of The Investigation 

As in our previous studies of EAS arrival directions, the main data sets consist of 513,602 showers 
registered with the EAS MSU array in 1984-1990 and 1,342,340 showers of the PRO-1000 array 
registered in 1997-1999. A description of the arrays can be found in [12] and [13] respectively. All 
EAS selected for the analysis satisfy a number of quality criteria and have zenith angles 9 < 45.7°. 
Showers registered with the two arrays have different number of charged particles N^ in a typical 
event. For the EAS MSU array, the median value of A'e is of the order of 1.6 x 10^, while that for 
the PRO-1000 array equals 3.7 x 10^. According to the modern models of hadronic interactions and 
data on the chemical composition of cosmic rays, the given values of N^ correspond to the energies 
of a primary proton Eq ~ 1.7 x 10^^ eV and Eq « 4.6 x 10^^ eV respectively, with an error of about 
10%-20%. The EAS MSU array allowed one to determine arrival directions of EAS with a better 
accuracy than the PRO-1000 array but the mean error is estimated to be of the order of 3° in both 
cases. 

The investigation is based on the method by Alexandreas et al., which was developed for the 
analysis of arrival directions of EAS registered with the CYGNUS array [14] and since then has been 
used multiple times for the analysis of results of other experiments on CRs. The idea of the method 
is as follows. To every shower in the experimental data set, arrival time of another shower is assigned 
in a pseudo-random way. After this, new equatorial coordinates (a, 5) are calculated for the "mixed" 
data set thus providing a "mixed" map of arrival directions. The mixed map differs from the original 
one but has the same distribution in declination 5. In order to compare both maps, one divides them 
into sufficiently small "basic" cells. A measure of difference between any two regions (cells) of the two 
maps located within the same boundaries is defined as 



S = (iVrcal - iVbg)/ViVb^, 

where A^rcai ^^nd A'^bg are the number of EAS inside a cell in the real and "mixed" (background) 
maps correspondingly. The mixing of the real map is performed multiple times, and the mixed maps 
are averaged then in order to reduce the dependence of the result on the choice of arrival times. 
The method is based on an assumption that the resulting mean "background" map has most of the 
properties of an isotropic background, and presents the distribution of arrival directions of cosmic rays 
that would be registered with the array in case there is no anisotropy. Thus, deviations of the real 
map from the background one may be assigned to a kind of anisotropy of arrival directions of EAS 
registered with the array. As a rule, selection of regions of excessive flux is performed basing on the 
condition 5 > 3. 

As in [15], the basic cells have the size 0.5° x 0.5° and the number of cycles of mixing equals 10,000. 
Due to this, the difference in the number of showers in basic cells of any two sequentially averaged 
maps by the end of the mixing cycle is of the order or (1-2) x 10^^ for both data sets. 

The procedure of constructing cells (regions) for comparing the number of EAS falling into them is 
as follows. Adjacent basic strips in 6 (each 0.5° wide) were joined into strips of width A6 = 3° . . . 30° 
with the step equal to 0.5° . Each wide strip was then divided into adjacent cells of some fixed width Aa. 
After this, we calculated the number of EAS inside each of these cells for both experimental (A'^rcai) 
and background (A^bg) rnaps. For each pair of cells, S was calculated then. A cell was considered as 
a possible cell of excessive flux (CEF) if S* > 3. Every wide strip was divided into cells with a shift of 
0.5° in the angular distance in a until all possible locations on the grid were covered. Contrary to [11], 



we did not restrict the search to cehs such that Aa = A(5/ cos S, where 6 is the mean value of 6 for the 
current strip, and Aa is rounded to the nearest half-integer number, but allowed Aa to deviate from 
that value in such a way that the area of the resulting cells deviated from the original one by 1/6. As 
a result, additional cells were considered. For example, besides cells of the size Aa x A6 = 3° x 3°, 
we also considered cells of the sizes 2.5° x 3° and 3.5° x 3°. Finally, only cells with more than 100 real 
EAS inside were analyzed. 

The method of Alexandreas et al. does not provide a direct answer to the question about the chance 
probability of appearance of a CEF. One is expected to calculate the corresponding probability basing 
on the value of S, which is an estimate of the standard deviation of a sample and thus acts as a 
significance level. It is implicitly assumed by this that the deviation of Aureal from Abg has a Gaussian 
distribution. Hence, the chance probability of a CEF to appear is estimated to be less than 1 — 0.9973 
providing it was selected at significance level S" > 3. 

In order to estimate the chance probability of appearance of a CEF basing on the number of 
EAS inside, we introduced the following simple method based on the binomial distribution [11]. Let a 
shower axis getting inside the CEF be a success. The number of trials equals the number of showers A^ 
in the data set under consideration, and an estimate of success (for a fixed region) equals p = Ny^g/N, 
where Abg is the number of showers in the cell of the background map. The assumption is based on 
the fact that in the method of Alexandreas et al., Abg is considered to be an expected number of 
showers in a cell. Obviously, the chance probability of finding exactly Afcai EAS in a cell (or region) 
equals 

where zv is a random variable equal to the number of successes in the binomial model, and C{N, Arcai) 
is the corresponding binomial coefficient. It is more interesting to consider a probability that there are 
at most Areai showers in the cell Pn = Pi^ < Areai)- The analysis performed and the data presented 
below demonstrate that values of Pn correlate well with the values of chance probabilities calculated 
on the basis of significance level S. 

The fact that S" > 3 for a given cell does not imply it is valid for its subcells. A situation such that 
S* < 3 for a number of adjacent cells but 5 > 3 for a cell (region) combined of them is possible. In 
order to improve the robustness of selection of CEFs, we have tried a number of additional quantities. 
One of the most efficient of them is the probability Pbc = P{i < -^bc)i calculated from the following 
binomial model. Let .^ be a random variable equal to the number of basic cells of the given cell with 
an excess of EAS over the background values, and N^^ is the corresponding experimental value. The 
number of trials equals the number of basic cells in the CEF. It is natural to assume the probability of 
success to be equal to 1/2. In the results presented below, all CEFs satisfy a condition Pbc > 0.9545, 
which corresponds to the significance level of 2(T for the Gaussian distribution. We thus reduced the 
chance that a CEF is selected solely due to the non-uniformity of the EAS distribution w.r.t. 5. 

3 Main Results and Discussion 

In what follows, we consider REFs found in the EAS MSU array experimental data under the conditions 

5>4, Pbc > 0.9545, Arcai>100. (1) 

Three regions satisfying (1) were found, see Fig. 1 and Table 1. Each of them consists of a number of 
partially overlapping cells, selected with the algorithm described above. The probability that any of 
the three REFs appears by chance is < 6 x 10~^. This estimate follows from the values of P/v shown 
in Table 1, as well as from the fact that S" > 4 for all the REFs: the corresponding probability equals 
K. 0.999937 for the significance level equal to 4. Regions A and B have "partners" in the PRO-1000 



data set. These are regions that satisfy conditions 5 > 3, Pbc 
overlap with those two regions in the EAS MSU data set. 
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Figure 1: Location of the regions of excessive flux in equatorial coordinates. Filled with red are the 
three regions found in the EAS MUS data that satisfy conditions (1). The open red rectangle shows 
the boundaries of a region found with the condition Pbc > 0.9545 omitted. Blue rectangles show 
REFs found in the PRO-1000 data set that partially overlap with the above REFs in the EAS MSU 
data and satisfy conditions S* > 3, Pbc > 0.9545, A^'rcai > 100 ("partners"). Filled red triangles show 
coordinates of the Galactic supernova remnants [16j. Filled blue diamonds show coordinates of the 
Galactic pulsars [17j . The curves at the left and right sides of the Figure show the Galactic plane, the 
n-like curve shows the Supergalactic plane. 



Table 1: Some parameters of the REFs in the EAS MSU data set that satisfy conditions (1) (A, B, C), 
and their "partners" in the PRO-1000 data set (A', B'). Notation: a and 6 are the ranges of values 
of the corresponding coordinates; A'rcah -^bg are the ranges of the experimental and background flux 
of EAS in the cells that constitute the REF respectively, max S and min P/v are the maximum value 
of the significance level S and the minimum value of Pn for the cells that constitute the REF. 



REF 


a, ° 


s,° 


A^rcal 


^bg 


max 5" 


minPjv 


A 
A' 


115.5. ..130.5 
112.5. ..136 


68... 72.5 
71. . . 78 


960. . . 1498 
2034. . . 6406 


841.9. . . 1350.8 
1881.4. ..6146.9 


4.233461 

3.852460 


0.999962 
0.999474 


B 
B' 


280. . . 287.5 
283. . . 292 


56. . . 60 

58... 62.5 


537. . . 1029 
1687. . . 3473 


451.5. ..904.3 

1552.4. ..3283.8 


4.478172 
3.556190 


0.999958 
0.999494 


C 


325... 329.5 


22. ..26.5 


108. . . 152 


67.34. . . 108.9 


4.954917 


0.999939 



The median values of A^e for the regions A, B, and C are close to the value obtained for the whole 
data set. The region that is found around region A when omitting the condition P^c > 0.9545 lies 
whithin the boundaries a = 96° . . . 130.5° and S = 67° . . . 79° and contains 7386 EAS. The median 



value of A'e for this REF is almost equal to the value for the whole data set. The minimum values 
of Pbc and P/v for the cells forming this region are equal to 0.831062 and 0.999962, respectively. 

As can be seen from the Figure, there are currently no known Galactic supernova remnants in the 
vicinity of the regions found. Some other possible accelerators of Galactic CRs, such as X-ray binaries 
and OB-associations have not been found there either. There is a bright open cluster NGC 2408 at 
the western boundary of Region A, but the distance from it to the Solar system is unknown, and it 
might be difficult to obtain a reliable estimation of its contribution to the flux of Galactic PeV CRs. 
On the other hand, there are isolated pulsars inside or in the close vicinities of all the regions. In 
view of this, the question arises as to whether these pulsars could accelerate charged particles to the 
energies of the order of PeV? 

To the best of our knowledge, there is currently no common model of acceleration of charged 
particles in the vicinity of isolated pulsars. To estimate the maximum energy of a particle accelerated 
near the light cylinder of a pulsar, we employ an expression that can be written as follows: -Emax = 
0.34 Z BQ,^ eV, where Z is the charge number of the particle, B is the strength of the surface magnetic 
field, G; Q is the angular velocity of a pulsar, rad s~^; and the radius of the pulsar equals 10^ cm [3]. 
The expression was obtained under the assumption that most of the magnetic field energy in the pulsar 
wind zone is transformed into the kinetic energy of the particles, and the density of electron-positron 
pairs does not exceed 10~^ of that of ions of iron. Table 2 presents approximate values of maximum 
energies of protons and iron nuclei accelerated by pulsars located in the vicinities of regions A, B, and 
C. It can be seen that in all three cases sufficiently heavy nuclei can be accelerated up to energies 
close to PeV. 

Table 2: Some parameters of the pulsars located near REFs A, B, and C, and estimates of the 
maximum energy of protons and iron nuclei. Notation: d is the distance from a pulsar to the Solar 
system, kpc; P is the barycentric period of rotation of a pulsar, s; B is the strength of the surface 
magnetic field, G; E is the rotation energy loss, ergs/s; maxE'p and maxE'pe are the maximum energies 
of a proton and an iron nucleus accelerated by a pulsar, eV. 



REF 


PSR Name 


d, 


P, 


Age, 


B, 


E, 


max Ep , 


max Epe , 






kpc 


s 


year 


G 


ergs/s 


eV 


eV 


A 


J0814+7429 


0.33-0.43 


1.292 


1.22 


10« 


4.72 


lO^i 


3.08 


103U 


4 


1012 


1 


1014 




J0700+6418 


0.48 


0.196 


4.52 


109 


1.17 


IQio 


3.61 


1030 


4 


1012 


1 


1014 




J0653+8051 


3.37 


1.214 


5.07 


W^ 


2.17 


1012 


8.37 


1031 


2 


1013 


5 


1014 




J0849+8028 


3.38 


1.602 


5.69 


10^ 


8.56 


10^1 


4.28 


1030 


4 


1012 


1 


1014 


B 


J1836+5925 


0.17-0.75 


0.173 


1.80 


10*^ 


5.20 


10^^ 


1.16 


10^4 


2 


10l4 


6 


lOl'^ 




J1840+5640 


1.70 


1.653 


1.75 


10^ 


1.59 


10^2 


1.31 


10^1 


8 


1012 


2 


1014 


C 


J2155+2813 


5.06 


1.609 


2.78 


10' 


1.23 


10^2 


8.68 


1030 


6 


10l2 


2 


10l4 




J2156+2618 


4.71 


0.498 


5.56 


10^ 


8.51 


IQio 


4.53 


1030 


5 


1012 


1 


1014 




J2151-F2315 


1.42 


0.594 


1.33 


10^ 


6.56 


lOii 


1.34 


1032 


2 


1013 


6 


1014 




J2139+2242 


4.71 


1.084 


1.21 


10^ 


1.26 


10^2 


4.41 


1031 


1 


1013 


4 


1014 



The distances to the pulsars from Table 2 exceed considerably the free path of PeV neutrons, 
and the amount of electromagnetic showers is insufficient for formation of REFs. In view of this, the 
question arises of how to keep the arrival direction of (charged) CRs propagating through the Galactic 
magnetic field. We are currently not ready to suggest a satisfactory model for answering this question, 
but we note at least two models were suggested to explain the REFs of CRs with energies of the order 
of 10 TeV registered with the Milagro gamma-ray observatory [18]. One of the models is based on CR 



acceleration in the vicinity of the Geminga pulsar [19], the properties of which are similar to those of 
the pulsar J1836+5925 located at the boundary of the B region [20]. The other one uses a specific 
magnetic field configuration able to focus a flux of charged particle [21]. 

To conclude, we have presented three regions of excessive flux of PeV CRs selected from the EAS 
MSU experimental data at a significance level > 4(T and satisfying several additional conditions. In 
the vicinity of every region, there are isolated Galactic pulsars that are able to accelerate sufficiently 
heavy nuclei to the required energies. In our opinion, this does not give enough ground to claim that 
the above REFs are formed due to contribution from these pulsars but can witness in favor of the 
assumption that (1) some isolated pulsars contribute to the flux of PeV CRs more than is usually 
assumed and (2) there are efficient mechanisms of focusing CR in the interstellar space. 

Only free, open source software was used for the investigation. In particular, all calculations were 
performed with GNU Octave [22] running in Linux. This research has made use of the SIMBAD 
database, operated at CDS, Strasbourg, France (http://simbad.u strasbg.fr/simbad). 

The research was partially supported by the Russian Foundation for Fundamental Research grant 
No. 08-02-00540 
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